Compare CFU and FungaLight for linear dilution of H2O2

require(tidyverse)
require(flowCore)
require(flowClust)
require(openCyto)
require(ggcyto)
require(cowplot)
#require(ggrdiges)

old <- theme_set(theme_minimal(base_size = 14) + panel_border(color = "gray20"))

Background & Goal

In the log dilution experiment, we compared FungaLight to the results of PI and CFU for C. glabrata cells treated with 0, 10, 100 and 1,000 mM hydrogen peroxide. We noted there that FungaLight seems to provide better separation in its score between the 100 mM and 1M treated sample than CFU and PI did. This led us to question whether FungaLight can be useful in distinguishing the levels of damage caused by different sublethal doses of hydrogen peroxide, for which CFU would give non-distinguishable readings.

Hanxi therefore performed an experiment on 2024-04-05 to 04-10, using a linear concentration series of hydrogen peroxide. This notebook records the pre-processing of that data and plotting the results.

C. glabrata post-treatment cells were stained with FungaLight (PI+SYTO9), and run through flow cytometry. The same sample was also plated for CFU.

Import data

FCS files are stored in rdss/Users/Hanxi/flow cytometry under the above dates. They are copied to this repository under the same names.

# use relative path to make it easier for collaboration
data.path = "../input/20240405-10 linear H2O2 correlation Cg"
dat0 <- read.flowSet(path = data.path, pattern = "*.fcs",
                     transformation = FALSE,  # the original values are already linearized. 
                     emptyValue = FALSE,  alter.names = TRUE,   # change parameter names to R format
                     column.pattern = ".H|FSC|SSC") # only load the height variables for the fluorescent parameters

Simplify the sample names

oriNames <- sampleNames(dat0)
tmp <- str_split(oriNames, pattern = "[ _]+", simplify = TRUE)[,c(1, 5)]
colnames(tmp) <- c("Date", "Treatment")
#treatment.levels <- c("mock", "5C", "10", "100", "1000")
sample <- data.frame(tmp) %>% 
  mutate(
    Date = mdy(Date),
    Treatment = gsub("\\.fcs", "", Treatment), # remove suffix
    Treatment = fct_inseq(Treatment),
    name = paste(format(Date, "%Y%m%d"), paste0(Treatment, "mM"), sep = "_")
  )
rownames(sample) <- oriNames
dat <- dat0 # make a copy
pData(dat) <- sample
write_tsv(pData(dat), file = "../input/20240405-linear-h2o2-cg-sample.tsv")

Gatting strategies

Load all gates from specification file

source("../script/20250206-FGL-flow-gate-specification.R")

Gating for singlets

The following gaphing steps are used to gate singlets by FSC and SSC values. Only singlets are included in analysis. Gate for outlier

#outlier.gate <- rectangleGate(filterId = "-outlier", "FSC.H" = c(1.2e5, 1e6), "SSC.H" = c(1e2, 1e6))
ggcyto(dat[1], aes(x = FSC.H, y = SSC.H), subset = "root") +
  geom_hex(bins = 64) + geom_gate(outlier.gate) + facet_wrap(~name, ncol = 2) + ggcyto_par_set(limits = "instrument")

Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Add gate to GS

gs <- GatingSet(dat) # create a GatingSet
# rename the samples
sampleNames(gs) <- pData(gs)$name
# add the outlier gate
gs_pop_add(gs, outlier.gate, parent = "root")

[1] 2

recompute(gs)

done!

Gate for singlets

scPars <- ggcyto_par_set(limits = list(x = c(0,1e6), y = c(30,300)))
ex <- Subset(dat[[1]], outlier.gate)
#polygon <- matrix(c(1e5, 1e5, 1e6, 1e6, 
#                    60, 105, 135,60), ncol = 2)
#colnames(polygon) <- c("FSC.H", "FSC.W")
#singlet.gate <- polygonGate(filterId = "singlet", boundaries = polygon)
ggcyto(ex, aes(x = FSC.H, y = FSC.W)) + geom_hex(bins = 128) + geom_gate(singlet.gate) + geom_stats() + scPars

Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Add this gate to the gatingSet

gs_pop_add(gs, singlet.gate, parent = "-outlier", name = "singlet")

[1] 3

recompute(gs)

done!

Gate for live cells

This is the live population with red fluorescence below 10^2.2

Main takeaway

1. Different from the last (log-dilution series) experiment, where we observed a large difference in %live between the mock and 100 mM treated cells, in this experiment the difference was not only variable across replicates but also much more modest, with the biggest change from 92% to 56%

2. %live drops to zero practically at 0.5M or above.

scPars <- ggcyto_par_set(limits = list(x = c(0,10^5.5), y = c(0,10^5.5)))
dates <- unique(pData(gs)$Date)
#polygon <- matrix(c(0, 5*10^3, 5*10^3, 0, # BL1.H, green
#                    10^2.2, 10^2.2, 0, 0),# BL3.H, red
#                  ncol = 2)
#colnames(polygon) <- c("BL1.H", "BL3.H")
#live.gate <- polygonGate(filterId = "live", .gate = polygon)
p.axis <- list(
  scale_x_logicle(breaks = 10^c(2,3,4,5)),
  scale_y_logicle(breaks = 10^c(2,3,4,5)),
  theme(
    strip.text = element_text(size = rel(1.1), face = 2)
  )
)

for(date in dates){
  p <- ggcyto(gs[pData(gs)$Date == date], aes(x = BL1.H, y = BL3.H), subset = "singlet") + 
    geom_hex(bins = 128) + 
    geom_gate(live.gate) +
    geom_stats(location = "data", adjust = c(0.005, 5), digits = 1) + 
    labs(title = paste0("Date: ", date)) +
    facet_wrap(~ fct_inseq(Treatment)) +
    p.axis + scPars
  print(p)
}

Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Add this gate to the gatingSet

gs_pop_add(gs, live.gate, parent = "singlet", name = "live")

[1] 4

recompute(gs)

done!

Gate for the intermediate population with high SYTO9

This is the intermediate population with high Red and Green. We hypothesize that they represent oxidatively damaged cells with partially compromised plasma membranes, resulting in more SYTO9 accumulation but not a significant increase in PI.

Main takeaway

1. For some reason, the percentage reported by geom_stats() is way off for 0.5 - 1M treated samples

2. % intermediate population peaked at 250 mM dose. Again, from 0.5M and up, we don’t see much distinction among the samples.

#polygon <- matrix(c(10^5, 1*10^3, 10^5,
#                    10^2.2, 10^2.2, 10^4), ncol = 2)
#polygon <- matrix(c(10^3, 10^5.5, 10^5.5,   # BL1.H, green
#                    10^2.2, 10^2.2, 10^4.5),# BL3.H, red
#                  ncol = 2)
#colnames(polygon) <- c("BL1.H", "BL3.H")
#inter.gate <- polygonGate(filterId = "inter", .gate = polygon)
for(date in dates){
  p <- ggcyto(gs[pData(gs)$Date == date], aes(x = BL1.H, y = BL3.H),
              subset = "singlet") + 
    geom_hex(bins = 128) + 
    geom_gate(inter.gate) + 
    geom_stats(location = "data", adjust = c(0.05, 5)) + 
    labs(title = paste0("Date: ", date)) +
    facet_wrap(~ fct_inseq(Treatment)) +
    p.axis + scPars
  print(p)
}

Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Plot to illustrate staining pattern

Show an example of the staining pattern with one replicate

scPars <- ggcyto_par_set(limits = list(x = c(0,10^5), y = c(0,10^4)))
p.ex <- ggcyto(gs[pData(gs)$Date == "2024-04-10"], 
               aes(x = BL1.H, y = BL3.H), subset = "singlet") + 
  geom_hex(bins = 64) + 
  geom_gate("live", colour = "gray20") + 
  #geom_stats(location = "data", adjust = c(0.05, 5)) + 
  labs(title = NULL) +
  facet_wrap(~ fct_inseq(Treatment)) +
  p.axis + scPars

Coordinate system already present. Adding new coordinate system, which will replace the existing one.

p.ex
ggsave("../output/20250120-linear-h2o2-example-FGL-staining.png", 
       plot = p.ex + labs(title = NULL, x = NULL, y = NULL) +
         theme(strip.text = element_text(size = rel(0.8), face = 1),
               axis.text = element_text(size = rel(0.8))),
       width = 3.6, height = 3.2)

Extract gated stats and MFI

% of events in each of the three gates

# we are only interested in the populations after singlet gating
nodes <- c("singlet", "intermediate", "live")
# get the event count
gated_cnt <- gs_pop_get_stats(gs, nodes, type = "count")
# get the MFI for each parameter in each population of interest
gated_mfi <- gs_pop_get_stats(gs, nodes, type = pop.MFI) %>% 
  select(sample, pop, starts_with("FungaLight"))

Combine the data and meta data

gated_stats <- full_join(gated_cnt, gated_mfi, by = c("sample", "pop")) %>% 
  right_join(pData(gs), by = c("sample" = "name")) %>% 
  relocate(Date, Treatment, .after = sample) %>% 
  relocate(sample, .after = last_col()) %>% 
  arrange(Date, as.numeric(Treatment))
write_tsv(gated_stats, file = "../output/20240405-linear-h2o2-Cg-gated-stats.tsv")

# we don't need the MFI for the latter analysis. add meta data to gated_cnt
gated_freq <- left_join(sample, gated_cnt, by = c("name" = "sample")) %>% 
  select(-name) %>% 
  group_by(Date, Treatment) %>% 
  mutate(perc = count / count[pop == "singlet"], .after = count,
         Treatment = fct_inseq(Treatment)) %>% 
  arrange(Date, Treatment)

#gated_stats <- gs_pop_get_stats(gs, type = "count") %>% 
#  as_tibble() %>% 
#  mutate(pop = gsub(".*/", "", pop), pop = gsub("-outlier", "cells", pop)) %>% 
#  pivot_wider(names_from = pop, names_prefix = "n_", values_from = count) %>% 
#  mutate(
#    #p_singlet = n_singlet / n_cells,
#    p_live = n_live / n_singlet,
#    p_intermediate = n_intermediate / n_singlet,
#    p_dead = 1 - p_live - p_intermediate
#  )
#write_tsv(gated_data, file = "../Input/20240124-fig-2-Dye-Order-gated-ht.tsv")
#print(gated_data)

final <- select(pData(dat), sample = name, everything()) %>% 
  left_join(select(gated_stats, -starts_with("n_")),  by = "sample") %>% 
  left_join(popMFI, by = "sample")

write_tsv(final, file = "../output/12.03-gated-output-20240214.tsv")
print(final %>% select(-sample))

Compare with CFU

Read CFU data

Here are the raw cfu data

cfu_raw <- read_tsv("../input/20240214-0405-h2o2-Cg-CFU-raw.tsv",
                    col_types = "cccccciiii", comment = "#")
cfu_calc <- cfu_raw %>% 
  select(Date, Treatment = H2O2_mM, Group, Dilution, Total) %>% 
  mutate(Date = ymd(Date),
         Treatment = fct_inseq(Treatment),
         CFU = Dilution * Total) %>% 
  # we only need the linear range experiment
  dplyr::filter(Date > ymd("20240401")) %>% 
  # group by Date to calculate the ratio of MO to MM
  group_by(Date) %>% 
  mutate(p_live = CFU / CFU[Group == "MM"])

Calculate survival from FungaLight

We calculate the percent live as the FGL survival. We also record the percent intermediate events and will investigate its utility for distinguishing doses not distinguishable by CFU.

dye_based <- gated_freq %>% 
  # record the total number of events as the count of singlets
  group_by(Date, Treatment) %>% 
  mutate(total = count[pop == "singlet"]) %>% 
  # now we don't need the singlet population and the count variable
  select(-count) %>% 
  dplyr::filter(pop != "singlet") %>%
  #ungroup() %>% 
  # we will pivot this data frame to a wider format, making 
  # two new variables out of the percentage variable
  pivot_wider(names_from = pop, values_from = perc, names_prefix = "p_") %>% 
  mutate(p_dead = 1 - p_live - p_intermediate) %>% 
  relocate(p_live, .after = total)

PI <- final %>% 
  dplyr::filter(Dye == "PI") %>% 
  select(Date, Treatment, Method = Dye, Score = p_live)

FGL <- final %>% 
  dplyr::filter(Dye == "Both") %>% 
  mutate(
    wFGL = p_live * 2/log(BL1_live, 10)
  ) %>% 
  select(Date, Treatment, wFGL, FGL = p_live) %>% 
  pivot_longer(cols = c(wFGL, FGL), names_to = "Method", values_to = "Score")
  

all <- bind_rows(
  PI, FGL,
  mutate(cfu, Method = "CFU") %>% rename(Score = CFU)
) 

Plotting the results

CFU survival

p.cfu <- cfu_calc %>% 
  ggplot(aes(x = Treatment, y = p_live)) +
  geom_bar(stat = "summary", fun = "mean", fill = "gray70", width = 0.8) +
  geom_point(aes(shape = as.character(Date)), 
             position = position_jitter(0.1),
             size = 1.5) + 
  scale_shape_manual(values = 4:6, guide = "none") +
  scale_y_continuous(labels = scales::percent, 
                     limits = c(NA, 1.05), breaks = seq(0,1,0.5)) +
  labs(title = "CFU", x = bquote(H[2]*O[2]~(mM)), y = "% Live") +
  theme(strip.text = element_text(size = rel(1.2), face = 2),
        strip.background = element_blank(),
        plot.title = element_text(hjust = 0.5))

p.cfu

FungaLight results

p.dye <- dye_based %>% 
  ggplot(aes(x = Treatment, y = p_live)) +
  geom_bar(stat = "summary", fun = "mean", fill = "gray70", width = 0.8) +
  geom_point(aes(shape = as.character(Date)), 
             position = position_jitter(0.1),
             size = 1.5) + 
  # overlay the p_intermediate points
  stat_summary(aes(y = p_intermediate), geom = "point", shape = 95,
            fun = "mean", size = 5, color = "steelblue3") +
  stat_summary(aes(y = p_intermediate, group = NA), geom = "line",
            fun = "mean", color = "steelblue") +
  geom_point(aes(y = p_intermediate, shape = as.character(Date)),
             position = position_jitter(0.1),
             color = "steelblue", size = 1.5, stroke = 1.5) + 
  scale_shape_manual(values = 4:6, guide = "none") +
  scale_y_continuous(labels = scales::percent, 
                     limits = c(NA, 1.05), breaks = seq(0,1,0.5)) +
  labs(title = "FungaLight", x = bquote(H[2]*O[2]~(mM)), y = "% Live") +
  theme(strip.text = element_text(size = rel(1.2), face = 2),
        strip.background = element_blank(),
        plot.title = element_text(hjust = 0.5))

p.dye

Alternative bar plots

p.dye.side_by_side <- dye_based %>% 
  select(Date, Treatment, p_live, p_intermediate) %>% 
  pivot_longer(cols = c(p_live, p_intermediate), 
               names_to = "par", values_to = "perc") %>% 
  mutate(par = factor(par, levels = c("p_live", "p_intermediate"),
                      labels = c("% live", "% intermediate"))) %>% 
  ggplot(aes(x = Treatment, y = perc, group = par)) +
  geom_bar(aes(fill = par), position = position_dodge(0.9),
           stat = "summary", fun = "mean", width = 0.8) +
  geom_point(aes(shape = as.character(Date), color = par), 
             position = position_jitterdodge(jitter.width = 0.1,
                                             dodge.width = 0.9),
             size = 1.2) + 
  scale_fill_manual(NULL, values = c("gray70", "lightblue2")) +
  scale_color_manual(NULL, values = c("black", "steelblue")) +
  scale_shape_manual(values = 4:6, guide = "none") +
  scale_y_continuous(labels = scales::percent, 
                     limits = c(NA, 1.05), breaks = seq(0,1,0.5)) +
  labs(title = "FungaLight", x = bquote(H[2]*O[2]~(mM)), y = NULL) +
  theme(strip.text = element_text(size = rel(1.2), face = 2),
        strip.background = element_blank(),
        plot.title = element_text(hjust = 0.5),
        legend.position = "inside",
        legend.position.inside = c(0.8, 0.9))

p.dye.side_by_side

Assemble the plot

plot_grid(
  p.cfu + labs(x = NULL, y = NULL) +
    theme(axis.text.x = element_blank(),
          title = element_text(size = rel(0.7)),
          axis.text = element_text(size = rel(0.7))),
  p.dye.side_by_side + labs(x = NULL, y = NULL) +
    theme(title = element_text(size = rel(0.7)),
          axis.text = element_text(size = rel(0.7)),
          legend.position = "none"),
  nrow = 2, rel_heights = c(1, 1.2)
)
ggsave("../output/20250120-linear-h2o2-CFU-FGL-compare.png",
       width = 3.1, height = 3.1)

Statistical tests

We use Tukey HSD test to compare all pairs of conditionsn for CFU

tmp <- aov(p_live ~ Treatment, data = cfu_calc)
TukeyHSD(tmp)

  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = p_live ~ Treatment, data = cfu_calc)

$Treatment
                  diff         lwr         upr     p adj
100-0    -9.352686e-01 -0.96550789 -0.90502935 0.0000000
250-0    -9.997259e-01 -1.02996516 -0.96948663 0.0000000
500-0    -9.999891e-01 -1.03022834 -0.96974980 0.0000000
750-0    -9.999891e-01 -1.03022834 -0.96974980 0.0000000
1000-0   -1.000000e+00 -1.03023927 -0.96976073 0.0000000
250-100  -6.445728e-02 -0.09469654 -0.03421801 0.0001301
500-100  -6.472045e-02 -0.09495972 -0.03448118 0.0001250
750-100  -6.472045e-02 -0.09495972 -0.03448118 0.0001250
1000-100 -6.473138e-02 -0.09497065 -0.03449211 0.0001248
500-250  -2.631738e-04 -0.03050244  0.02997609 1.0000000
750-250  -2.631738e-04 -0.03050244  0.02997609 1.0000000
1000-250 -2.741035e-04 -0.03051337  0.02996516 1.0000000
750-500   9.251971e-17 -0.03023927  0.03023927 1.0000000
1000-500 -1.092966e-05 -0.03025020  0.03022834 1.0000000
1000-750 -1.092966e-05 -0.03025020  0.03022834 1.0000000

For p_live with FungaLight

tmp <- aov(p_live ~ Treatment, data = dye_based)
TukeyHSD(tmp)

  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = p_live ~ Treatment, data = dye_based)

$Treatment
                  diff        lwr         upr     p adj
100-0    -0.1298926649 -0.3394148  0.07962944 0.3562860
250-0    -0.8678501009 -1.0773722 -0.65832800 0.0000001
500-0    -0.9469346756 -1.1564568 -0.73741257 0.0000000
750-0    -0.9478779476 -1.1574001 -0.73835584 0.0000000
1000-0   -0.9481361662 -1.1576583 -0.73861406 0.0000000
250-100  -0.7379574360 -0.9474795 -0.52843533 0.0000007
500-100  -0.8170420107 -1.0265641 -0.60751991 0.0000002
750-100  -0.8179852826 -1.0275074 -0.60846318 0.0000002
1000-100 -0.8182435013 -1.0277656 -0.60872140 0.0000002
500-250  -0.0790845747 -0.2886067  0.13043753 0.7961753
750-250  -0.0800278466 -0.2895500  0.12949426 0.7885837
1000-250 -0.0802860653 -0.2898082  0.12923604 0.7864874
750-500  -0.0009432719 -0.2104654  0.20857883 1.0000000
1000-500 -0.0012014906 -0.2107236  0.20832061 1.0000000
1000-750 -0.0002582187 -0.2097803  0.20926389 1.0000000

For p_intermediate with FungaLight

tmp <- aov(p_intermediate ~ Treatment, data = dye_based)
TukeyHSD(tmp)

  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = p_intermediate ~ Treatment, data = dye_based)

$Treatment
                  diff         lwr         upr     p adj
100-0     0.1124892688 -0.08745318  0.31243172 0.4516826
250-0     0.5711061447  0.37116370  0.77104859 0.0000066
500-0    -0.0415257799 -0.24146823  0.15841667 0.9787052
750-0    -0.0450275307 -0.24496998  0.15491492 0.9699612
1000-0   -0.0452750554 -0.24521750  0.15466739 0.9692633
250-100   0.4586168759  0.25867443  0.65855933 0.0000631
500-100  -0.1540150487 -0.35395750  0.04592740 0.1739532
750-100  -0.1575167995 -0.35745925  0.04242565 0.1589407
1000-100 -0.1577643242 -0.35770677  0.04217813 0.1579233
500-250  -0.6126319246 -0.81257437 -0.41268948 0.0000031
750-250  -0.6161336754 -0.81607612 -0.41619123 0.0000029
1000-250 -0.6163812001 -0.81632365 -0.41643875 0.0000029
750-500  -0.0035017508 -0.20344420  0.19644070 0.9999999
1000-500 -0.0037492755 -0.20369172  0.19619317 0.9999998
1000-750 -0.0002475247 -0.20018997  0.19969492 1.0000000

plot_grid(
  p.cfu + scale_x_discrete(limits = rev) + coord_flip(),
  p.dye + scale_x_discrete(limits = rev) + coord_flip() + 
    theme(axis.title.y = element_blank(), axis.text.y = element_blank()),
  rel
)

---

all %>% 
  dplyr::filter(!Treatment %in% c("mock", "5C"), Method != "wFGL") %>% 
  mutate(Method = fct_relevel(Method, "CFU", "FGL", "PI")) %>% 
  ggplot(aes(x = Treatment, y = Score)) +
  geom_bar(stat = "summary", fun = "mean", fill = "gray") +
  #geom_point(size = 2, position = position_jitter(0.05)) + 
  stat_summary(fun.data = "mean_sdl", fun.args = list(mult = 1), 
               geom = "errorbar",  width = 0.1) +
  labs(x = bquote(H[2]*O[2]~(mM)), y = "% Live") +
  facet_wrap(~Method, nrow = 1) +
  theme_cowplot() +
  theme(strip.text = element_text(size = rel(1.2), face = 2),
        strip.background = element_blank())
ggsave("../output/20240323-for-R35-CFU-compared-with-FGL-and-PI.png", width = 5, height = 3)

Example of sample-to-sample variability, using 100 mM treated sample as case

p <- ggcyto(gs[pData(gs)$Treatment %in% c("100", "250")],
            aes(x = BL1.H, y = BL3.H), subset = "singlet") +
  geom_hex(bins = 64) +
  geom_gate(live.gate) +
  geom_stats(location = "gate", adjust = c(2, 0.5), digits = 1) +
  facet_grid(Treatment ~ Date) + 
  labs(title = NULL) +
  p.axis + scPars

Coordinate system already present. Adding new coordinate system, which will replace the existing one.

p
ggsave("../output/20250214-example-variability-between-replicates.png", 
       width = 5, height = 3.5)

---

Below are hanxi’s plotting code

combined_data <- read_tsv("/space/htang5/Documents/FUN1_data_analysis/biol-4386-course-project-htang5/Methods Paper Writing/02-FungaLight-Distinguishing-Power/Input/20240219-PI-Both-Distinguishing-Power-gated-CFU-ht.tsv")
print(combined_data)
combined_data %>% 
  dplyr::filter(Dye %in% c("Both")) %>% 
  dplyr::filter(Treatment %in% c("0", "10", "100")) %>% 
  ggplot(aes(x = FGS, y = CFU_Survival)) +
  geom_point(aes(color = Treatment))

combined_data %>% 
  dplyr::filter(Dye %in% c("Both")) %>% 
  ggplot(aes(x = Treatment, y = FGS)) +
  geom_point(aes(color = Date))

combined_data %>% 
  dplyr::filter(Dye %in% c("Both")) %>% 
  ggplot(aes(x = Treatment, y = percent_inter)) +
  geom_point(aes(color = Date))
#ggsave("../output/20240104-PI-Dilutions-Gated_RP-Plot.png", width = 7, height = 5)

combined_data$new_FGS <- 1*combined_data$percent_Live - 0.6*combined_data$percent_inter - combined_data$percent_Dead
pf1g <- combined_data %>%
  pivot_longer(cols = c(CFU_Survival, new_FGS),
               names_to = "Assay",
               values_to = "Survival") %>% 
  dplyr::filter(!(Dye %in% c("PI"))) %>% 
  dplyr::filter(!(Treatment %in% c("5C")))  %>%
     mutate(Assay = factor(Assay, levels = c("CFU_Survival", "new_FGS"), labels = c("CFU", "FGL"))) %>% 
  #mutate(Dilution = factor(Dilution, levels = c("250 x", "500 x", "1000 x", "2000 x"))) %>%
  ggplot(aes(x = Assay, y = Survival)) +
  facet_wrap(~ Treatment, nrow = 1, labeller = as_labeller(c("0" = "0 mM", "10" = "10 mM", "100" = "100 mM", "1000" = "1000 mM"))) +
  stat_summary(fun.data = "mean_cl_boot", geom = "pointrange", color = "red2",
               size = 0.8, position = position_nudge(x = 0.2)) +
  geom_point(stroke = 1, size = 2, position = position_jitter(0.01)) +
  scale_y_continuous() +
  theme_minimal()

pf1g +
  labs(x = "Assay Type", y = "Survival Quantification") +
  panel_border(color = "black", size = 1.5) +
  theme(axis.line = element_blank(),
        strip.background = element_blank(),
        axis.text.x = element_text(size = 20, face = "bold", color = "black"),  # Adjust size and face (bold)
        axis.text.y = element_text(size = 20, face = "bold", color = "black"),
        axis.title.x = element_text(size = 20, face = "bold"),
        axis.title.y = element_text(size = 20, face = "bold")
        )+ 
  theme(strip.text = element_text(size = 20, face = "bold")) 
ggsave("../output/20240301-Cg-Log-Gated_.png", width = 7, height = 5)

combined_data$percent_PIDead <- 1 - (combined_data$n_PIDead / combined_data$n_singlet)
p <- combined_data %>%
  pivot_longer(cols = c(CFU_Survival, percent_PIDead),
               names_to = "PIAssay",
               values_to = "PISurvival") %>% 
  dplyr::filter(!(Dye %in% c("Both"))) %>% 
  dplyr::filter(!(Treatment %in% c("5C")))  %>%
     mutate(PIAssay = factor(PIAssay, levels = c("CFU_Survival", "percent_PIDead"), labels = c("CFU", "PI"))) %>% 
  #mutate(Dilution = factor(Dilution, levels = c("250 x", "500 x", "1000 x", "2000 x"))) %>%
  ggplot(aes(x = PIAssay, y = PISurvival)) +
  facet_wrap(~ Treatment, nrow = 1, labeller = as_labeller(c("0" = "0 mM", "10" = "10 mM", "100" = "100 mM", "1000" = "1000 mM"))) +
  stat_summary(fun.data = "mean_cl_boot", geom = "pointrange", color = "red2",
               size = 0.8, position = position_nudge(x = 0.2)) +
  geom_point(stroke = 1, size = 2, position = position_jitter(0.01)) +
  scale_y_continuous() +
  theme_minimal()

p +
  labs(x = "Assay Type", y = "Survival Quantification") +
  panel_border(color = "black", size = 1.5) +
  theme(axis.line = element_blank(),
        strip.background = element_blank(),
        axis.text.x = element_text(size = 20, face = "bold", color = "black"),  # Adjust size and face (bold)
        axis.text.y = element_text(size = 20, face = "bold", color = "black"),
        axis.title.x = element_text(size = 20, face = "bold"),
        axis.title.y = element_text(size = 20, face = "bold")
        )+ 
  theme(strip.text = element_text(size = 20, face = "bold")) 
ggsave("../output/20240301-Cg-Log-Gated-PI.png", width = 7, height = 5)

pf1g <- combined_data %>%
  dplyr::filter(!(Dye %in% c("Both"))) %>% 
  dplyr::filter(!(Treatment %in% c("1000", "5C")))  %>%
  #mutate(Dilution = factor(Dilution, levels = c("250 x", "500 x", "1000 x", "2000 x"))) %>%
  ggplot(aes(x = Treatment, y = CFU_Survival)) +
  facet_wrap(~ Dye, scale = "free_y") +
  geom_boxplot() +
  geom_point(stroke = 1, size = 2, position = position_jitter(0.2)) +
  scale_y_continuous(labels = scales::percent) +
  theme_minimal()
pf1g +
  labs(x = "Hydrogen Peroxide Treatment Concentration", y = "CFU Survival") +
  panel_border(color = "black", size = 1.5) +
  theme(axis.line = element_blank(),
        strip.background = element_blank(),
        axis.text.x = element_text(size = 12, face = "bold", color = "black"),  # Adjust size and face (bold)
        axis.text.y = element_text(size = 16, face = "bold", color = "black"),
        axis.title.x = element_text(size = 16, face = "bold"),
        axis.title.y = element_text(size = 16, face = "bold")
        )+ 
  theme(strip.text = element_text(size = 16, face = "bold")) 
#ggsave("../output/20240104-PI-Dilutions-Gated_RP-Plot.png", width = 7, height = 5)

—>